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ABSTRACT 

The cold dark matter (CDM) structure formation model predicts that about 5-10 % of a typical 
galactic halo of mass ~ 1O 12 M 0 is in substructures with masses < 1O 8 M 0 . To directly detect such 
substructures, we propose to observe dust continuum emission from a strongly lensed QSO-host galaxy 
using a large submillimeter interferometer array with a high angular resolution of ~ O.Olarcsec such as 
the planned Atacama Large Submillimeter Array (ALMA). To assess their observational feasibility, we 
numerically simulate millilensing of an extended circular source by a CDM substructure modeled as 
a tidally truncated singular isothermal sphere (SIS) embedded in a typical QSO-galaxy lens system, 
B1422+231, modeled as a singular isothermal ellipsoid (SIE) with an external constant shear and 
a constant convergence. Assuming an angular resolution of O.Olarcsec, we find that the angular 
positions of ~ 1O 8 M 0 substructures at several kpc from the center of the macrolens halo can be 
directly measured if the size of the dust continuum emission region and the gradient of the surface 
brightness at the position of the perturber are sufficiently large. From the astrometric shift on a scale 
of a few times 10 mas of an image perturbed by a subhalo with respect to an unperturbed macrolensed 
image, we can break the degeneracy between subhalo mass and distance provided that macrolensing 
parameters are determined from positions and fluxes of multiple images. 

Subject headings: cosmology: theory - dark matter - gravitational lensing - large-scale structure of 
universe 


1. INTRODUCTION 

Recent high-resolution iV-body simulations of the cold 
dark matter (CDM) structure formation model show an 
excess in the number of subhalos (“substructures”) in 
a galaxy-sized halo in comparison with the number of 
Milky Way or M31 satellites (Klypin et al. 1999; Moore 
et al. 1999). If the CDM-based structure formation sce¬ 
nario is the correct one, then we expect a large number 
of invisible satellites inside a halo of galaxies. 

Currently, gravitational lensing is the only probe that 
can test the abundance of these invisible satellites or 
equivalently, substructures in a galaxy-sized halo. In 
fact, gravitationally multiply lensed QSOs have recently 
been used for putting limits on the surface density and 
the mass of such invisible substructures (Mao & Schnei¬ 
der 1998; Metcalf & Madau 2001; Chiba 2002; Metcalf & 
Zhao 2002; Dalai & Kochanek 2002; Bradac et al. 2002). 
Statistical analyses on these quadruple QSO-galaxy lens¬ 
ing systems that show an anomalous flux ratio imply the 
presence of substructures along the line-of-sight to the 
multiple images (Metcalf & Zhao 2002; Evans & Witt 
2003; Keeton, et al. 2003). 

However, the location of a subhalo along the line-of- 
sight is poorly constrained. In fact, more massive ex- 
tragalactic 10 8 — 1O 9 M 0 halos along the line-of-sight to 
the multiple images may also significantly contribute to 
altering the flux ratio of multiple images (Metcalf 2005). 
In order to determine the distances to these subhalos, we 
need to have more information in addition to the position 
and flux of multiple images (Yonehara, et al. 2003). 

To directly map out the three-dimensional locations 
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of 1O 8 M 0 invisible subhalos, we propose to observe the 
dust continuum emission from a strongly lensed QSO- 
host galaxy using a large submillimeter interferometer 
array with a high angular resolution of ~ 0.01 arcsec 
such as the ALMA. 

Recent submillimeter observations revealed that about 
60 % of strongly lensed active galactic nuclei (AGNs) are 
luminous at a submillimeter wavelength of 850/rrn (Bar- 
vainis & Ivison 2002) implying that a large amount of en¬ 
ergy emitted by stars and AGNs in the QSO host galaxy 
has been absorbed by dust grains and then re-radiated 
at longer wavelengths. The dust is heated to tempera¬ 
tures of 20-50K and radiates as a modified blackbody at 
far-infrared wavelengths (Wiklind 2003). At a redshift 
0 — 3 — 4, the corresponding energy peak is at a submil¬ 
limeter wavelength ~ 0.1 mm. The region responsible for 
this dust continuum emission has a typical scale of 10 1 2 pc 
to a few kpc. Therefore, the effects of an extended source 
with an angular size larger than the Einstein radius of a 
perturber will be easily observed at submillimeter wave¬ 
lengths (Blain 1999; Wiklind & Alloin 2002). 

In this paper, we show our simulation results of 
millilensing of an extended circular source by a CDM 
substructure whose Einstein radius is smaller than the 
source size. We also show that the astrometric shifts 
of images perturbed by a subhalo can reveal the three- 
dimensional position of a sublralo along the line-of-sight. 
As an example, we choose a typical QSO-galaxy lensing 
system, B1422+231, with cusp-caustics in order to check 
the observational feasibility of mapping substructures. 

2. SIMULATION 

The typical QSO-galaxy lensing system B1422+231 
shows an anomalous flux ratio. The images consist 
of three highly magnified images A, B, and C, and a 
faint one D located near the lens galaxy. The observed 
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mid-infrared flux ratios are A/\B\ = 0.94 ± 0.05 and 
C/\B\ = 0.57 ± 0.06 (Chiba et al. 2005). The ra¬ 
dio and optical flux ratios are consistent with the mid- 
infrared values within a ~ 10 % error. If we assume that 
the gravitational potential of the macrolens is a smooth 
one, we expect a cusp-caustic relation in the flux ra¬ 
tio, / = (A + B + C)/{\A\ + |B| + |Cj) = 0. However, 
B1422+231 shows a 20 % deviation from the / = 0 re¬ 
lation. Assuming that such an anomaly is caused by the 
presence of an SIS perturber in the macrolens halo, the 
anomaly might be in image A (Dobler & Keeton 2005). 
The redshifts of the source and the lens are zs = 3.62 and 
zl = 0.34, respectively. In what follows, we assume the 
following cosmological parameters: the present density of 
matter fl m = 0.3, the present density of the cosmological 
constant = 0.7, and the Hubble parameter h = 0.7, 
which yield the angular diameter distances to the lens 
and to the source, Dl = 1.00 Gpc and Ds = 1.49 Gpc, 
respectively. 

To model the macrolens system, we adopt an SIE in 
an external shear field in which the isopotential curves 
in the projected surface perpendicular to the line-of-sight 
are ellipses (Kormann et al. 1994). For further details, 
see Chiba (2002). In what follows, we assume an angular 
resolution of 0.01 arcsec which will be achieved by the 
planned ALMA at 850 /tm. Because the Einstein radius 
for an SIS with a one-dimensional velocity dispersion a 
at Dl = l.OOGpc and Ds = 1.49Gpc is 

Qev~ 1-0 x lof-mas, ( 1 ) 

P \21 km s -1 J 


observation with an angular resolution of O.Olarcsec will 
reveal subhalos with a one-dimensional velocity disper¬ 
sion (7 > 20 km s -1 for B1422+231. 

As a model of a CDM subhalo, we consider a 
spherically symmetric, tidally truncated SIS with one¬ 
dimensional velocity dispersion er. At a distance r from 
the center of a macrolens galactic halo with a one¬ 
dimensional velocity dispersion cto, the tidal radius is 
approximately given by r t « ra/ao. For an SIS with a 
one-dimensional velocity dispersion a ~ 20 km s” 1 , the 
tidally truncated mass is Mgjg ~ 10 8 Mg. Provided that 
the mass function for subhalos satisfies dn/dM oc M -2 , 
as many iV-body simulations suggest (e.g., Klypin et al. 
1999), the differential lensing cross section dS per log¬ 
arithmic mass interval is dS = 9 E dn/dlnM oc M 1 / 3 


because the area inside an Einstein radius for a tidally 
cut SIS is 6 2 Ep oc Mgjg. Thus, the contribution of mas¬ 
sive subhalos to the substructure lensing is somewhat 
significant in comparison with that of less massive ones. 
We also assume that a perturbing subhalo is located at 
a distance equal to the size of the Einstein radius of the 
macrolens halo ve = DlOe = 3.8kpc from the macrolens 
center (where 0e is the Einstein angular radius of the 
macrolens halo). This assumption can be verified as fol¬ 
lows. Consider a column with a unit section centered 
at a point P aligned to the line-of-sight direction cen¬ 
tered at a distance equal to the size of the Einstein ra¬ 
dius of the macrolens re from the macrolens center. We 
introduce a column distance L from P as a coordinate 
in the direction of the line-of-sight in the column. If 
the subhalo distance is larger than the size of the Ein¬ 
stein radius of the macrolens halo, r rg, the tidally 


truncated mass approximately satisfies Mgjg oc L , since 
r ~ L. Therefore, the probability of having a subhalo 
with a fixed one-dimensional velocity dispersion a per 
logarithmic column distance InL is dp(L) ~ L~ 1 dlnL, 
provided that the mass function for subhalos in the col¬ 
umn with a fixed one-dimensional velocity dispersion sat¬ 
isfies dn/dM oc M~ 2 . Thus, the lensing probability is 
expected to be larger for a smaller column distance L. 

If the “real” mass function takes a different form, or the 
above assumption is incorrect, the perturbing subhalo 
may reside at a position with a larger column distance 
L > te- Then the corresponding subhalo mass will be 
increased by a factor of ~ L/te- 

In our simulation, we put a tidally truncated SIS with 
a mass of 2 x 1O 8 M 0 near the A image . 3 The corre¬ 
sponding Einstein radius of the SIS is 9e p ~ 17 mas. 
As a model of emission from cold dusts, we assume a 
circularly symmetric Gaussian source with a standard 
deviation L = 2.5 x 10 2 pc. 

3. RESULT 

The left panels in Figure 1 show lensed images with¬ 
out any perturbers, whereas the right panels show those 
perturbed by an SIS perturber. As shown in a zoomed in 
picture of image A ( top right), the angular position of an 
SIS can be clearly identified by a dipole structure that 
consists of a pair of dark and bright spots at the place 
where the surface brightness gradient is non-vanishing 
(Inoue & Chiba 2003). One can also notice that shifts 
of an image outside the Einstein radius 9e p of an SIS 
with respect to the unperturbed macrolensed image are 
significantly suppressed ( top right). 

If we can approximate the lens system near the per¬ 
turber locally as an SIS plus a constant shear 7 and a con¬ 
stant convergence k, then the lens equation normalized 
by the Einstein radius 9e p of an SIS in the coordinates 
aligned to the shear can be written in terms of a source 
position y = ( 2/1 , J/ 2 ) and an image position x = (aq,^) 
as 

y = (l-r)x-A (2) 

where 

r =('A 1 A)- ( 3 > 

Let Xo = (xoi, X 02 ) be the image position for a macrolens 
model without any perturbers. In terms of the polar 
coordinates (R p , </>) defined as x = (R p cos <f>, R p sin </>) , 
the lens equation © yields 

*oi(^ P -^" 1 ) 2 +*o2(^-^r 1 ) 2 = 

( 4 ) 

where £1 = 1 — k — 7 and £2 = 1 — k + 7 . The astro- 
metric shift Ax = x — xo as a function of xo can be 
obtained by solving the fourth-order equation 0 with 
equation m analytically. For a positive parity case £1 > 
0 , £2 > 0 , an image at a horizontal axis (* 01 , 0 ) is shifted 
to (a:oi + Cf 1 ,!)), whereas an image at a vertical axis 

3 We choose an appropriate position of the SIS that approxi- 
mately gives the observed anomaly flux ratio, but we did not fix 
the relative position of the images in a small source size limit. 
Therefore, our simulation should be interpreted as a demonstra¬ 
tion of extended source size effects rather than that of a best-fit 
model. 
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Fig. 1.— Simulated images of B1422+231 at submillimeter wavelengths without any perturbers (left) and with an SIS perturber (right) 
for a Gaussian circular source with a standard deviation L = 2.5 x 10 2 pc. A tidally truncated SIS with a mass of 2 x 1O 8 M0 is put at the 
center of a thick circle (top right) near the A image. The dashed circle represents a circle with the tidal radius centered at an SIS perturber. 
The angular resolution is assumed to be 0.01 arcsec. Astrometric shifts in the tangential direction and those in the radial direction in the 
coordinates aligned to the shear are represented by arrows. 


(0, £ 02 ) is shifted to (0, £02 + £^’ 1 ) (Inoue & Chiba 2004). 
One can see in the top right panel in figure 1 that im¬ 
age A, with positive parity k ~ 0.38,7 ~ 0.47, is largely 
stretched in the tangential direction by ^-j -1 ~ 6.7 and 
slightly stretched in the radial direction by £2 ~ 0.92 

with respect to the center of the perturber. The astro- 
metric shift is significantly suppressed if the lensed im¬ 
age is placed outside the tidal radius centered at the SIS 
perturber, because the projected gravitational potential 
decreases faster, as R~ l , than an SIS without trunca¬ 
tion, where R is the projected distance from the center 
of the SIS. If the redshift and the one-dimensional veloc¬ 
ity dispersion uo of the macrolens halo are known, then 
we obtain the tidal radius as a function of the distance 
from the center of the macrolens halo. From the observed 
angular size for which the astrometric shifts are signif¬ 
icantly suppressed, we can determine the distance from 
the center of the lens halo, in principle. 


Next, we show that a measurement of astrometric 
shifts can break the degeneracy between substructure 
mass and distance in the line-of-sight to the image. Con¬ 
sider a halo at redshift z = Zh and a clump (a substruc¬ 
ture) in the foreground of the halo at redshift z c < Zh- 
Let D s , D c h, D cs , and Dh be the angular diameter dis¬ 
tances to the source, between the clump and halo, be¬ 
tween the clump and source, and to the halo, respec¬ 
tively. Provided that the angle between the perturber 
and the macrolensed image is sufficiently small in com¬ 
parison with the Einstein radius of the macrolens 9e, the 
effective lens equation for an SIS plus a constant shear 
and a constant convergence can be written as (Keeton 
2003) 

y=(l-f)x-^, (5) 

where y = (1 —/3r) _1 y, f = 1 — (1 —/3r) _1 (l — T), and 



4 


INOUE and CHIBA 


(3 = D c hD s /D cs Dh- Because an effective unperturbed 
image position x 0 that satisfies y = (1 — r)x 0 is equal 
to the unperturbed macrolensed image position x 0 , as¬ 
trometric shifts can be written as Ax = x — Xq. The 
effective convergence and shear in T are then (Keeton 
2003) 

_ (i -/3 )(k -/?(k 2 - q 2 )) _ (1 — /3 )t 

KeS (1 - /3k) 2 - 0V ’ 7eff (1 - /3k) 2 - /3V' 

(6) 

Thus, astrometric shifts are described in the same man¬ 
ner as in the case when (3 = 0. A tangential shift 
is Axi = (1 — Kgjj — 7 e ff) , and a radial shift is 
Ax 2 = (1 — K e ff + 7 e ff ) -1 (Inoue & Chiba 2004). As 
shown in figure 2 , for images with positive parity, the 
tangential shift decreases as the distance between a per- 
turber and a macrolens halo increases. As the distance 
to the clump D c decreases, the astrometric shift becomes 
less conspicuous and more isotropic. In fact, in the limit 
D c —* 0 or, equivalently, [3 —> 1, the astrometric shifts 
converge to a perfect isotropic shift Ax = x 0 /|x 0 |. Thus, 
the distance between a halo and a clump can be mea¬ 
sured by the absolute value and the anisotropy of the 
shifts, provided that the distances to the macrolens halo 
Dh and the source D s are known already. 

Similarly, we can calculate astrometric shifts Ax in 
the background case z c > Zh■ Again, the effective lens 
equation takes the same form as in equation © (Keeto n 
2003), 

y=(l-f)x-|p (7) 

where x = (1 - /3T)x, f = 1 - (1 T)(l - /3T)~\ 
and (3 = D c hD s /D c D s h ( D s h is the angular diame¬ 
ter distance between the source and the halo). We 
can define an effective position of the unperturbed im¬ 
age Xo satisfying y = (1 — r)x 0 , which leads to Xq = 
(1 — /3T)x 0 . Then astrometric shifts can be expressed as 
Ax = (1 — /3r) - 1 Ax((l — /3r)x 0 ), where Ax = x — x 0 . 
Thus, astrometric shifts in the background case are dif¬ 
ferent from those in the foreground case with the same 
/ 3 , although the effective lens equation takes the same 
form in both cases. From the right panels in Figure 2 
and 3, one can see that the shifts of images near the hor¬ 
izontal axis (xi, 0 ) are larger for the background case in 
comparison with the foreground case with the same j3. 
As shown in the left panel in Figure 2 and both panels 
in Figure 3, the astrometric shifts of images on the hori¬ 
zontal axis Ax = (£)“ , 0) and those on the vertical axis 
Ax = (O,^ 1 ) do not depend on f3 in the background 
case. However, astrometric shifts of all the other images 
do depend on (3. As (3 increases from 0, astrometric shifts 
tend to become isotropic except for images on the axes 
of coordinates aligned to the shear, but they cannot be 
perfectly isotropic as in the foreground case even if f3 = 1. 

Now let us assess our assumption for constancy in the 
convergence k and shear 7 . The absolute errors in as¬ 
trometric shifts in the horizontal direction (aligned to 
the shear) and in the vertical direction are written in 
terms of absolute errors in the convergence 8k and shear 
Sj as <5(Axi) ~ '?r 1 (l^ K / K l + \&l/l\) and 8(Ax 2 ) ~ 
^2 i\^ K / K \ + I 67 / 7 I), respectively. In order to achieve 
a three-dimensional mapping with an angular resolution 


of a unit Einstein radius (~ 9e p = 1), the errors 8k and 

8j should satisfy I4 _1 |(I ( W K I + |^ 7 / 7 |) l(i = 1 , 2 ) 
at the scale of a unit Einstein radius. In order to check 
this condition, we compute k and 7 in the neighborhood 
of an SIS perturber for our SIE model. As shown in 
Figure 4, the variances in k and 7 at the scale of the 
perturber’s Einstein radius 9e p = 17 mas are only of the 
order of a few %. Thus, our approximation is verified for 
an SIS with one-dimensional dispersion <j ~ 20 km s -1 . 
In more general settings, the magnitude of spatial vari¬ 
ation in the convergence and shear at the scale of the 
Einstein radius 9e p of a perturber can be estimated as 
Sk/k ~ < 57/7 ~ 9e p /9e■ If the macrolens is also an SIS 
with a one-dimensional velocity dispersion a 0 , then we 
should have ^(cr /a o ) 2 C 1. In our model, this condition 
is roughly satisfied, because cro ~ 180 km s _1 . However, 
for subhalos with a much larger velocity dispersion er, we 
should take the effect of spatial variance in the conver¬ 
gence k and shear 7 into account. 

4. OBSERVATIONAL FEASIBILITY 

The observed values of the flux density at 850pm of 
QSOs at 2 = 3 — 4 with a starburst are typically several 
tens of mjy, which implies an emission from dust with a 
temperature of 30-60K from a fairly extended region of 
10 2 — 10 3 pc (Barvainis & Ivison 2002). Let us assume 
that the far-infrared dust emission can be represented 
by a single-temperature graybody emission with a rest- 
frame spectral energy distribution (SED) (Wiklind 2003) 

S(u r ,T d ) = Q,B(v r ,T d )(l - exp(-(iz r /r' 0 ) /3 )), ( 8 ) 

where B(v r , T d ) is the Planck function in terms of a rest- 
frame frequency v r and a dust temperature T d , (3 = 1 — 2 
is the model parameter that controls the grain emissivity, 
uq is the critical frequency that controls the opacity of 
the dust, and 0 is the solid angle of the source seen by an 
observer. A rest-frame SED S is related to the observer- 
frame SED as S(v) = S((l + z)v, T d )/(1 + z). Assuming 
a model with (3 = 1.5, v$= 6 THz, a dust temperature 
T d =40 K, and a dust emission region from a circular 
disk with a radius L = 250 pc, at 850 /im, the observed 
energy density would be ~ 2.6 mjy. 

Because the ratio of the area D within the Einstein 
radius of a 2 x 10 s M 0 SIS perturber to the area of im¬ 
age A is approximately 1/50, the observed energy flux 
from the area D will be ~ 4 x 10 -4 Jy, taking a mag¬ 
nification factor of ~6 into account. In order to mea¬ 
sure the difference in the energy flux at the scale of the 
Einstein radius of the SIS perturber, for a source with 
a smooth surface brightness distribution (e.g., Gaus¬ 
sian) one needs a sensitivity of at least ~ 4 x 10~ 5 Jy. 
On the other hand, the sensitivity of ALMA at 850 
pm for a point source (signal-to-noise ratio S/N= 1) is 
1.2 x 10 -5 x (integration time/5hour) - 1 / 2 Jy. 4 Thus, 
ALMA has a sufficient sensitivity for a source with a 
spatially varying surface brightness to detect the dipole 
structure caused by a 2 x 1O 8 M 0 SIS perturber. The re¬ 
quired integration time would be several hours. Intrinsic 
surface brightness variations in the source might make 
the identification of the clear dipole structure rather dif¬ 
ficult. However, such a difficulty can be avoided by 

4 See the ALMA home page at 

http: / / www.eso.org/projects/alma/science/bin / sensitivity.html 
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comparing multiple macrolensed images, because intrin¬ 
sic variations in one macrolensed image can be linearly 
mapped to those in another one. 

5. SUMMARY AND DISCUSSION 

In this paper, we have proposed a new method to real¬ 
ize three-dimensional mapping of CDM substructures in 
extragalactic halos. We have shown that at submillime¬ 
ter wavelengths, a measurement of astrometric shifts of 
perturbed multiple images with respect to unperturbed 
images can break the degeneracy between subhalo mass 
and position in the line-of-sight to the image if resolved at 
the scale of an Einstein radius of the perturber. Further¬ 
more, assuming that the tidal radius is approximately 
equal to the subhalo size, the distance from the center 
of the macrolens halo can be determined from the size of 
the region within an image outside which the astrometric 
shift is suppressed. 

Although we have considered a very simple model, the 
essential features of astrometric shifts will not be dra¬ 
matically altered even if we consider a more complicated 
model such as a macrolens perturbed by a number of 
non-circular subhalos and/or subhalos having a so-called 
NFW profile (Navarro, Frenk & White 1996). This is 
because (1) the astrometric shifts caused by a signifi¬ 
cant number of smaller mass subhalos are canceled out 
on larger scales if they are not clustered on that scale. 
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However, several massive perturbers (often not within 
the image) may distort the images so that the accuracy 
in measuring the distance to a single perturber along 
the line-of-sight may become significantly worse. (2) De¬ 
viation from a circular symmetry or from an SIS mass 
profile can affect the astrometric shifts significantly. For 
instance, if the substructures had an NFW or some other 
less-cuspy inner profile, the simulated dipole pattern in¬ 
side the Einstein radius would be less prominent (Inoue 
& Chiba 2004). However, the generic feature of sub¬ 
structure lensing, i.e., its dependence on the parity of a 
macrolensed image outside the Einstein radius of a per¬ 
turber, will probably be unaltered; if the outer density 
profile of a substructure falls off as fast as an SIS with in¬ 
creasing radius, then we expect similar astrometric shifts 
to those for an SIS. Even if the outer density profile is 
shallower than an SIS, we can make a distinction be¬ 
tween models with different density profiles and ellip- 
ticities from astrometric shifts of a number of internal 
positions within the source, if the source has a substruc¬ 
ture. 
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Fig. 2.— Astrometric shifts Ax for an SIS perturber in a macrolens halo for (3 = 0 {left) and those for an SIS perturber in the foreground 
of a macrolens for f3 = 0.5 (right) with respect to images for an unperturbed macrolens in the image plane (x\,X 2 ) aligned to the shear. 
Both horizontal and vertical axes are normalized by the Einstein radius 6e p of an SIS perturber. We adopt a convergence k = 0.38 and a 
shear 7 = 0.47 for the macrolens. 



Fig. 3.— Astrometric shifts Ax for an SIS perturber in the background of a macrolens halo for (3 = 0.9 (left) and (3 = 0.5 (right) with 
respect to images for an unperturbed macrolens model. The other parameters are the same as those in Figure 2. 
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Fig. 4.— Contour plots of k (top) and 7 (bottom) for our SIE model of B1422+231 in the neighborhood of an SIS perturber (the position 
is plotted in the top right panel in Figure 1) that is centered at the coordinate center. 






